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Abstract 

We consider a model proposed earlier by us for describing a form of plastic insta- 
bility found in creep experiments . The model consists of three types of dislocations 
and some transformations between them. The model is known to reproduce a number 
of experimentally observed features. The mechanism for the phenomenon has been 
shown to be Hopf bifurcation with respect to physically relevant drive parameters. 
Here, we present a mathematical analysis of adiabatically eliminating the fast mode 
and obtaining a Ginzburg-Landau equation for the slow modes associated with the 
steps on creep curve. The transition to the instability region is found to be one of 
subcritical bifurcation over most of the interval of one of the parameters while su- 
percritical bifurcation is found in a narrow mid-range of the parameter. This result 
is consistent with experiments. The dependence of the amplitude and the period of 
strain jumps on stress and temperature derived from the Ginzburg-Landau equation 
are also consistent with experiments. On the basis of detailed numerical solution via 
power series expansion, we show that high order nonlinearities control a large portion 
of the subcritical domain. 
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1 Introduction 



Instabilities in plastic flow has been an object of attention for a long time in metallurgical 
literature. Conceptually simplest form of this [Lubahn & Felgar, 1961; Hall, 1970] manifests 
when the material is subjected to a creep test wherein a force is applied and the response 
in the form of elongation of the material is measured, which under normal conditions, is a 
smooth strain-time curve. However, under certain metallurgical conditions, one sees steps 
on creep curve suggesting a form of instability [Lubahn & Felgar, 1961; Hall, 1970; Da 
Silveira & Monteiro, 1979]. A better known form of the instability having the same physical 
origin, but conceptually more difficult, appears when the material is deformed under tensile 
deformation [Bodner & Rosen, 1967; Brindley & Worthington, 1970; Penning, 1972]. Here, 
the material is subjected to a predetermined response namely a constant rate of deformation, 
and the force or the stress developed in the sample is sought to be measured. Even in this 
case, one finds a smooth stress-strain curve under normal conditions. However, when the 
system is in the regime of instability, for some values of material parameters the stress-strain 
curve exhibits multiple load drops. Each of the load drops is related to the formation and 
propagation of dislocation bands [Chihab et al, 1987]. It is in the latter type of testing 
where plastic instability manifests much more easily than in the former. The phenomenon is 
referred to as the Portevin-Le Chatelier (PLC) effect or the jerky flow and is seen in several 
metals such as commercial aluminium, brass, alloys of aluminium and magnesium [Brindley 
& Worthington, 1970]. It is observed only in a window of strain rates and temperature. 
In contrast, the phenomenon of steps on creep curve, which is the subject of the present 
discussion, is seen in few instances [Da Silveira & Monteiro, 1979; Zagarukuyko et al, 1977; 
Stejskalova et al , 1981 ]. The reason attributed to this is that it is difficult to obtain proper 
control on metallurgical parameters wherein this form of instability can be observed. There 
is one more form of instability which manifests under constant stress test where similar 
strain jumps are seen as a function of time. The origin of instability in all these three modes 
of testing is known to be the same.( For a current status of both experiments and theory 
see Kubin et al, 1993). Eventhough steps on creep curve are seen in a limited number of 
experiments, it is often straight forward to translate many of the experimental results of 
constant strain rate test into that of creep test. It is generally agreed that the microscopic 
origin of the instabilities arises due to the interaction of dislocations with mobile point defects 
and is referred to as dynamic strain ageing. This leads to negative strain rate characteristic. 
The basic idea was formulated by Cottrell [1953] few decades ago. 

The early phenomenological models are all static since they do not deal with time de- 
velopment. In contrast, methodology of dynamical systems addresses precisely this aspect. 
Until about a decade and a half ago, there were no models which looked at the problem 
from the point of view of bifurcation theory. An attempt to understand the problem in the 
above perspective was first made by our group several years ago [Ananthakrishna & Sahoo, 
1981b; Ananthakrishna & Valsakumar, 1982, 1983; Valsakumar & Ananthakrishna, 1983]. 
In a series of papers starting from an extended Fokker-Planck equation for the distribution 
function of the velocity of dislocation segments, and then splitting this into a mobile and an 
immobile component, we arrived at a model which consisted of three types of dislocations 
and some transformations between them [Ananthakrishna & Sahoo, 1981a, 1981b; Sahoo 
& Ananthakrishna, 1982; Valsakumar & Ananthakrishna, 1983]. The basic idea could be 
summarized by stating that the phenomena is due to a Hopf bifurcation resulting from 
nonlinear interactions between three different types of dislocations, suggesting a new math- 
ematical mechanism for the instability. Eventhough the spatial inhomogeneous structure 
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was ignored and only the temporal oscillatory state was sought to be described, the model 
and its extensions to the case of constant strain rate test, proved to be very successful in 
that it could explain most of the experimentally observed features (for critical reviews see 
Kubin & Lepinoux, 1988 and Kubin et al, 1993). For instance, the existence of negative 
strain rate sensitivity in a regime of strain rates, comes out naturally as a consequence of 
the Hopf bifurcation [Ananthakrishna & Valsakumar, 1982]. It must be emphasized, that 
this feature has been measured in most experiments on the PLC effect, and is assumed in all 
other theoretical models [ Cottrell, 1953; Bodner & Rosen, 1967; Penning, 1972; Kubin and 
Estrin 1990; see also Kubin & Martin, 1988; Martin & Kubin, 1992; Kubin et al 1993]. Other 
experimentally observed features such as the existence of bounds on strain rate for the PLC 
effect to occur, the existence of critical strain and its dependence on applied strain rate, the 
dependence of the amplitude on the strain etc., also follow. One other important prediction 
which is direct consequence of the dynamical basis of the model is the existence of chaotic 
stress drops in a window of strain rates [Ananthakrishna & Valsakumar, 1983; Ananthakr- 
ishna & John, 1990]. Recently this prediction has been verified which in turn implies that 
only a few degrees of freedom are required for a dynamical description of the phenomenon. 
This offers justification for the use of only a few degrees of freedom for the description of the 
temporal aspect, eventhough the system is spatially extended. [Ananthakrishna et at, 1995; 
Ananthakrishna & Noronha, 1995; Quaouire & Fressengeas, 1995; Venkadesan et al 1995; 
Noronha et al 1996 & 1997]. (Note that a spatially extended system implies infinite degrees 
of freedom.) Description of the phenomenon which includes the initiation and propagation 
of the bands during the PLC effect has also been recently attempted [Ananthakrishna, 1993]. 

Since the introduction of bifurcation theory into this field several years ago by our group, 
there has been a resurgence of interest in plastic instabilities in the light of introduction of 
new methodology borrowed from the theory of dynamical systems. This has further helped 
to obtain new insights hitherto not possible [Kubin & Martin, 1988; Estrin & Kubin, 1989; 
Kubin & Estrin, 1990; see also Kubin & Martin, 1988; Martin & Kubin, 1992; Kubin, 1993; 
Ananthakrishna et al, 1995a &b and references therein]. One of the aims of such theories is 
to be able to relate the microscopic dislocation mechanisms to the measurable macroscopic 
quantities. However, in the process, we feel that finer aspects of dynamical systems have 
been glossed over in this field. For instance, one often finds that casual remarks are made 
about fast and slow modes without actually going through the procedure of demonstrat- 
ing the existence of such modes and eliminating the fast modes in favour of the slow ones 
[Aifantis, 1988; Hahner & Kubin, 1992; Hahner, 1993]. In addition, under the adiabatic elim- 
ination, the slow modes are complicated functions of the original modes. Yet, hand waiving 
arguments have been used in building models which we believe are technically suspect. 

The purpose of this article is at least three fold. First, we wish to illustrate the power 
and utility of dynamical methods to the study of dynamical aspects of the model. We 
demonstrate the technical aspects of adiabatic elimination (largely addressed to metallurgists 
working in the area of plastic instabilities) and then derive the equation for the slow modes 
having the form of time-dependent Ginzburg-Landau (TDGL) equation. Second, we wish 
to compare the results of the dependence of the amplitude and the period of the strain 
jumps on stress and temperature with experiments on creep curve as also with the results 
translated from the constant strain rate case to the creep case. This will allow us to relate the 
theoretically introduced parameters to the macroscopically measured quantities. Third, we 
wish to compare these results and the detailed numerical results via power series expansion 
with the approximate solutions for the limit cycles obtained earlier [Ananthakrishna & Sahoo, 
1981b; Valsakumar & Ananthakrishna, 1983]. 
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The principal results of this work are as follows. We show that the derived quintic TDGL 
equation for the slow modes is valid in the mid-range of the instability domain. Within this 
range, there is a narrow sub-range of the parameter where the system exhibits a supercritical 
bifurcation, beyond which subcritical bifurcation is exhibited. The results of this calculations 
on the amplitude and period of the strain jumps are compared with experiments. In addition, 
analytical work based on the TDGL equation and a detailed comparative numerical study 
via the power series solution shows that the model exhibits unusual properties in a certain 
domain of the parameter. The study further shows that very high order nonlinearities 
actually govern the nature of the subcritical bifurcation beyond the domain handled by the 
quintic TDGL equation. 

The plan of the paper is as follows. In what follows (section 2) we present a brief summary 
of the model. In section 3, the technique of extracting the TDGL equation using reductive 
perturbative approach is applied to the model. In section 4, approximate limit cycle solution 
obtained through the TDGL equation is compared with experiments and with the numerical 
solution of the model as well as earlier studies on the model. Section 5 contains summary 
and discussion. 



2 A Model for Steps on Creep Curve 

We start with a brief summary of the model. The basic idea of the model is that instabilities 
in plastic flow are a consequence of nonlinear interactions between different types of dislo- 
cation populations. Spatial dependence is ignored and only temporal aspects are described, 
the idea being that the few degrees of freedom used for describing the phenomenon (disloca- 
tion populations) correspond to the collective modes of dislocations in the spatially extended 
system. As mentioned in the introduction, this finds support from the recent experimen- 
tal verification of chaos in the PLC effect. The details of the model can be found in the 
original references [Ananthakrishna & Sahoo, 1981b; Valsakumar & Ananthakrishna, 1983]. 
The model consists of mobile dislocations m and immobile dislocations im and another type 
which mimics the Cottrell's type i, which are dislocations with clouds of solute atoms. Let 
the corresponding densities be N m , N im and N i: respectively. The basic dislocation mecha- 
nisms included are the following: (a) production of dislocations by cross glide (m m + m 
, 6 is the cross glide coefficient, V m is the velocity of mobile dislocations; V m = V (cr a /a ) m 
[Alexander, 1986] where a a is the applied stress), (b) immobilization of two mobile disloca- 
tions (m + m — > im + im, (3 is the rate at which two mobile dislocations get immobilized) 
and annihilation of two mobile dislocations (m + m — > 0), and annihilation of a mobile with 

an immobile dislocation [m + im 0), (c) once a mobile dislocation is immobilized, it 
can be reactivated by athermal or thermal means (i m), (d) lastly, a mobile dislocation 
can acquire a cloud of solute atoms and move with them. We consider such dislocations as 
distinct from both the mobile as well as the immobile and denote them i. This process is 
represented by m ^5 i, a m being a function of the concentration of solute atoms. However, 
as the solute atoms gather progressively around dislocations, these dislocations eventually 
will be immobilized {i ^4 im). ( Note that in principle, f3, j3' and (3 are different.) These 
mechanisms lead to the rate equations for the densities of dislocations: 

N m = 9V m N m -(P + P')N^-0N m N im + -fN im -a m N m , (1) 
N im = PN^-pN im N rn -'yN im + a i N i , (2) 
Ni = a m N m - oijiVj. (3) 
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(For the sake of simplicity, we take f3' — and f3 — f3.) One can give a more transparent 
interpretation for Eq. (3) which represents whole process of acquiring a cloud of solute 
atoms by a mobile dislocation and eventual immobilization. This can be seen by defining 
Ni = f^K^t — t')N m (t')dt'i where the kernel K(t) which represents the fact that solute 
atoms are arriving at a certain rate has been chosen to be a m e~ ait . It can be easily checked 
that Eq. (3) is the differential form of the definition of iVj. One last remark about some of 
the terms appearing in these equations. The first two terms plus the last term in Eq.(l), and 
equivalent terms in Eqs.(2), were derived starting from a velocity distribution function of 
the dislocation segments [ Ananthakrishna & Sahoo 1981a; Sahoo & Ananthakrishna, 1982]. 
Thus, this model is an improved model to explain the instability in plastic flow. Equations 
(1-3) can be cast into a dimensionless form by using scaled variables: 



x = N m (~), V = N im (Jr),z = iVi(— ), t = 9V m t, (4) 
7 0V m 7« m 



to get 



x = (1 — a)x — b x 2 — xy + y, (5) 
y = b (b x 2 -xy-y + azj, (6) 
z = c(x — z), (7) 

where a = ^jtS fro = mr) an( ic = The dot represents differentiation with respect to r. 
Equations (5-7) are coupled set of nonlinear equations which support limit cycle solutions for 
a range of parameters a, bo and c, that are physically relevant, a refers to the concentration 
of the solute atoms, bo refers to the reactivation of immobile dislocations and c to the time 
scales over which the slowing down occurs. We have demonstrated the existence of limit 
cycle solutions and also obtained approximate closed form solutions for the limit cycles 
[Ananthakrishna & Sahoo, 1981b; Valsakumar & Ananthakrishna, 1983]. In addition, the 
model has been studied numerically. Using the Orowan equation which relates the rate of 
change of strain(S') to mobile dislocation density and the mean velocity: S = bN m V m , with 
b as the Burger's vector, steps on the creep curve follow automatically since the densities of 
dislocations are oscillatory. Several experimental results are reproduced [Ananthakrishna & 
Sahoo, 1981b; Valsakumar & Ananthakrishna, 1983]. 



3 Reductive Perturbative Approach 

For the sake of completeness, we outline the reductive perturbative approach. Near the point 
of Hopf bifurcation of the system (Eqs.5-7), corresponding to a critical value of the drive 
parameter, a pair of complex conjugate eigenvalues and another real negative eigenvalue 
exist for the linearized system of equations around the steady state. As we approach the 
critical value from the stable side, the real parts of the pair of complex conjugate eigenvalues 
approach zero and hence the corresponding eigenvectors get slower and slower. In contrast, 
the effect of the change in the drive parameter on the real negative eigenvalue is negligible. 
Thus, the eigenvector corresponding to the real negative eigenvalue is the fast mode. For this 
reason, the slow modes are the ones which determine the formation of new states of order. 
The reductive perturbative method is a method where the asymptotic equation is extracted in 
a systematic way [Taniuti & Wei, 1968; Newell & Whitehead, 1969; Kuramoto & Tsuzuki, 
1974; Mashiyama et al, 1975; Richter et al, 1981]. The method involves in first finding 
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the eigenvectors corresponding to the slow modes and looking for a solution in the subspace 
spanned by these vectors. The effect of the nonlinearity is handled progressively as a per- 
turbation of the linear solution in power series of the deviation from the critical value of the 
drive parameter. It may be worth emphasizing that this method is essentially the same as 
reduction to center manifold. Indeed, the equivalence of the center manifold theory [Carr, 
1981; Guckenheimer & Holmes, 1983; Troger & Steindl, 1991] with the reductive pertur- 
bation has been established [Chen et al, 1996]. Other techniques of extracting amplitude 
equations have been devised whose end results are basically the same. For instance, pertur- 
bative renormalization group method [Goldenfeld et al 1989; Chen et al, 1996] and its recent 
extension on the basis of envelope theory [Kunihiro, 1995 & 1996] has also been developed as 
a tool for global asymptotic analysis which can be used to extract the amplitude equations. 

In the present model, there is one simplification arising from the fact that x can be 
identified as the fast (stable) mode, and hence this will be eliminated adiabatically. Then 
we derive the equation for the slow dynamics through a power series expansion. Although 
this adiabatic elimination is approximate, the result is found to be in a reasonably good 
agreement with the accurate method (which will be reported elsewhere) except for finer 
details. 

Consider the Eqs.(5-7). It can be shown that there is a domain of instability for the 
parameters a, b and c. The range of values of b is ~ — 1CT 2 , that of c is ~ — 1CT 1 
while a ~ 1. Thus, b and c < 1. Thus, there are three time scales that are in principle 
different from each other. Physically, they correspond to the mobile (fastest a), the Cottrell 
type (c) and the immobile (slowest b ). Under these conditions we can eliminate one full 
equation adiabatically. This aspect becomes transparent if the above equations are written 
in terms of variables which are deviations from the steady state. There is only one steady 
state defined by: 

= = 1 - 2a+l(1 - ft 2a)2 + 86 ° |1/2 . and, = 1/2. (8) 

Defining new variables which are deviations from the steady state 

X = x - x a , Y = y - y a , and Z = z - z a , (9) 

Eqns. (5-7) take the form 

X = -(aX + X Y + b X 2 + XY), (10) 
Y = -b (TX + 5Y -aZ - b X 2 + XY) , (11) 
Z = c(X-Z). (12) 

where 

a = a + 2b x a + y a - 1, x = x a - 1, 



(13) 



r = y a - 2b x a ,5 = x a + 1. 
Now, we can rescale the time-like variable by r' = b r and get 



b o^7 = -(aX + X Y + b X 2 + XY), (14) 
^7 = - (VX + 5Y -aZ- b X 2 + XY) , (15) 
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It is clear that in the limit of small bo, \X\ — > oo, unless the left hand side vanishes identically. 
Thus, we can eliminate X in favour of the other two variables and obtain: 

dY 

— = -b [aX + 2(x a + X)Y - aZ] , (17) 
ar 

and — = c\X-Z], (18) 
ar 

where 

X = W (" (C * + Y) + ^ + 2(a " 2boX)Y + " 2] ") • (19) 

(The other root for X is unphysical since it corresponds to negative dislocation density.) 
It must be emphasized that this adiabatic elimination becomes more exact as the value of 
the parameter 6 gets smaller. Equations (17-18) will be solved reductive perturbatively. 
Writing these equations as a matrix equation where the nonlinear part appears separately 
from the linear part, we get 

^=LR + N (20) 
or 

where 

( V )> ( 21 ) 




(22) 

with Co = b (^ — 2x a ) and the nonlinear part, N, is given by 

# _ / -ab (X + IY) - 2b XY \ _ ( E^ =2 t n Y n \ ( . 

\ <X + *Y) )-\Y£*CnY»)- {Z6) 

The coefficients ^ n and ( n appearing in the last expression for iV are functions of a, b and 
c whose expressions for the required first few are given in Eqs.(A8-Al5) of Appendix A. 

Consider stability of the fixed point as a function of the parameter c. The eigenvalues, 
A±, of the matrix L are 

A± = - (cq - c ± yf((co - cf - 8b x a c)^j . (24) 

The fixed point becomes unstable when c is less than Co (c is non-negative) and the discrim- 
inant of Eq.(24) is negative giving a pair of complex conjugate eigenvalues. This holds when 
the inequality (3a — 1)(1 — 2a 2 ) > 2&o(2 + a) 2 is satisfied. The instability region in the - 
versus a plane for a fixed value of &o(— 10~ 4 ) is shown in Fig. 1. It extends approximately 
between a — | and a = 4j. 

To get approximate analytical solution of Eq.(20), we follow the reductive perturbative 
approach similar to that used by Mashiyama et al. [1975] and Richter et al. [1981]. We 
choose c = c (l — e) with < f C 1 and write the matrix L as a sum of two matrices, 
L = L + eLi, where L is the matrix L evaluated for c = c and 

(25) 
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The eigenvalues of L are A ± = where u = J (2b c x a ). Considering that the solution 
for R grows continuously out of the critical eigenmodes, we can express it as 

R( T ) = ^(r)e^ r r+ + **(r)e _,W7 "rl, (26) 

where \l/ and are time- dependent complex amplitudes while r± are right eigenvectors 
defined by L r± = ±zur±. Similarly, we introduce left eigenvectors, defined by s^L = 
±zu;s±. (Note that s^r_ = s^r + = 0. T refers to transpose operation of the matrix 
concerned). Substituting this expression for R in Eq.(^Of) and multiplying it by from the 
left-hand side gives the equation governing 

e^ T — = e^e lu,T - e^*e^ T + V PJ^ T + c.c.) n . (27) 

dT n=2 

fi — nx + tfi2 and P n are complex coefficients which are functions of a, &o and c. Their 
expressions are, respectively, given in Eqs. (Al) and (A4) in Appendix A. We express \I> in 
power series of 

* = e'^i + eip 2 + ••• (28) 

and introduce slowly varying multiple time scales t±(= er), r 2 (= e 2 r),... replacing j- by 
+ e-£- + e 2 -£^ + .... With this substitution in Eq. (|2T|) , it can be solved successively by 

equating terms of equal power in e on both sides of the equation. At first, O(e^) terms give 
rise to the equation 

w = (29) 

which means that, on the time scale of r, ipi is constant. The next higher order, 0(e), terms 
give the equation 

P 2 (^ie luJT + c.c.)V^ T (30) 



dr 

which can be integrated to get ip 2 - 0(e^) terms give the equation 



<9-?/>3 dtp 



+ ^ = //Vi - Mte- 2 *" + 2P 2 [(ij 1 ij 2 e 2lLUT + Vi^ 2 *) + cc] + P^ x e luJT + c.c.) 



dr <9ti 

(31) 

from which we extract the slow dynamics: 

+ ^1 + V\^i (32) 

Since to O(e^) ^ = e^ipi, Eq. (|32|) gives the dynamics for the complex amplitude 

d^ , , 9 , s 

— = e/i* + 7?|^| 2 $. (33) 
ar 

The expression for the complex coefficient rj = 771 + 2772 is given in Eq. (A2) in Appendix 
A. Equation (^) is a cu&ic TDGL equation where ^ is the complex amplitude for the slow 
mode which is also referred to as complex order parameter (of the new state of temporal 
order). Its steady state solution gives the amplitude squared as 

|H/| 2 = - e ^I (34) 



and the associated frequency, Q, (with \I> = |^/|e* nT ) as 

Sl = e[ii2- —Hi I • (35) 



Note that \l/ describes the limit cycle oscillation. This solution exists provided r]i is negative 
since ii\ is positive (see expression of fi in Eq.(Al) in Appendix A), rji is found to be 
negative in a narrow region (0.4713 < a < 0.5120) as can be seen in Fig. 2. In this case, 
since the amplitude of the slow modes grows continuously in proportion to (see Eq . ([34]) ) 
the transition is supercritical bifurcation (continuous or 'second order') For other values 
of a (within the instability range), r]i is positive implying that the transition is subcritical 
bifurcation (discontinuous or 'first order'). One has to then go to quintic or even higher terms 
in the TDGL equation for obtaining an expression for the slow modes (order parameter). 
We have carried out the reductive perturbative method further to derive the quintic TDGL 
equation: 

= e/i^ + r]\^\ 2 ^ + z/|^| 4 ^. (36) 

(XT 

The steady state solution of this equation gives the amplitude squared as 



*| 2 = ^ -- + J(-) 2 -4^ , (37) 




= ( V2 - Vi-) l^l 2 + e (to - to-) • (3* 



and the frequency as 

We found V\ to be negative in the range 0.4538 < a < 0.5097 (see Fig. 2). This means a 
narrow region: 0.4538 < a < 0.4713 within the subcritical bifurcation exists over which the 
quintic TDGL equation has a solution (see Fig. 2). 



4 Comparison with Experiments and Numerical Solu- 
tions 

4.1 Comparison with Experiments 

To start with, consider the supercritical regime where expressions for the amplitude and 
period of the limit cycle are simple. Using the steady state solution of the cubic TDGL 
equation, Eqs.(|3"4|) & (|3^), the dependence of |\l/| 2 and period, P ~ , of the limit cycle on 
a are, to the leading order, found to be a~ 2 and (const. +a 2 ), respectively. The parameter a(= 
a m /6V m ) is a function of the applied stress (a a ) and temperature (T). As remarked earlier, ot m 
is proportional to the concentration of the solute atoms, and therefore a m ~ exp(—E/kT). 
Using the standard expression [ Alexander 1986] V m (a a ,T) = Vo(cr a /cr ) m exp(— E m /kT) 
(with m > 1), we get 

1 / a„\ m E- Em 

e k T . (39) 



From this we see that the amplitude of the limit cycle has a strong increasing dependence on 
stress and an increasing dependence on temperature if we assume that E > E m . In contrast, 
the period of the limit cycle has a weak decreasing dependence on stress (as the leading 
contribution is constant) and a decreasing dependence on temperature. Since stress and 



9 



temperature are measurable quantities, our predictions can be compared with experimental 
results. The amplitude and period of the limit cycle are related, respectively, to the amount 
of strain jumps and the period of the jumps on the creep curve through the Orowan equation. 
There are very few experiments in this mode of testing as mentioned in the introduction. The 
only experiment where this dependence on stress and temperature has been measured is that 
by Zagoruyko et al [1977]. However, it is possible to translate the results from experiments 
in constant strain rate case to the creep case and compare them with the results of this 
calculation. According to Zagoruyko et al [1977], the amplitude of the strain jumps increases 
with stress while its period has a weak but decreasing dependence on stress. Experiments 
from constant strain rate case also exhibit the same trend when the results are translated in 
terms of constant stress experiments. It is well known that the amplitude of the stress drops 
decreases with applied strain rate. In fact, even the present model predicts this behaviour 
for the constant strain rate mode [Ananthakrishna & Valsakumar 1982]. This implies that 
the dependence of the amplitude of strain jumps on stress should increase [ Kubin et al 
1993]. [ This relation can be seen as follows. In constant strain rate case, the deformation 
rate is fixed and the stress developed in the sample is measured. When the contribution to 
the plastic strain rate increases due to increased dislocation motion (for whatever reasons), 
the stress has to fall in order to keep the applied strain rate constant. Thus, the relation 
between strain rate and stress is opposite.] Clearly, the general trend is consistent with the 
experimental results. Fig. 3 shows the actual dependence of the amplitude and the period 
on a in the domain of supercritical bifurcation (not the leading one as given by the above 
expressions). Zagarukuyko etal [1977] also report that the amplitude of the strain jumps 
increases while its period decreases with temperature, which is consistent with our result ( 
provided E > E m ). 

In the case of subcritical bifurcation, the dependence of the amplitude of the strain jumps 
and the period on a have similar behaviour as in the case of the supercritical bifurcation. 
Eventhough, in this case, leading order dependence on a is difficult to obtain we have their 
actual dependence plotted in Fig. 4. It may be noted that they have qualitatively the same 
behaviour as in Fig. 3 with the degree of dependence more pronounced in this case. 

Experiments in constant strain rate case show that the stress drops are seen to arise both 
abruptly as well as continuously [Kubin et al, 1988]. Translating this result to the constant 
stress case, it implies that the strain jumps can arise both abruptly and continuously. This 
feature again follows from our calculations. 

4.2 Comparison with Numerical Solutions 

Having derived the TDGL equation, we would first like to compare its result with the 
numerical solutions obtained via Eqs. (17-18) and via Eq. (|20|) . Secondly, due to the fact 
that the quintic TDGL equation (Eq. |36[) is valid in a limited domain, one suspects that 
higher order nonlinearities are controlling the subcritical bifurcation, which is quite unusual. 
Therefore, we would like to investigate the numerical solution obtained by solving Eq.(|2T)|) 
keeping successive leading powers in Y. This analysis confirms the above suspicion that 
higher order nonlinearities are important in this model outside the domain of validity of the 
quintic TDGL equation. 

Using the steady state solution of the quintic TDGL equation, Eq. (|36|) , in Eq.(p6l), we 
get an approximate expression for the limit cycle which is usually called the secular motion. 
The equations governing the secular motion are given by Eqs. (B5) and (B6) (hereafter 
called secular equations) in Appendix B. This can be compared with the numerical solutions 
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obtained using Eqs. (17-18) (hereafter called the reduced model) and with that obtained 
using Eq . (|20"D where the expansion in Y of Eq. (23) is limited up to a certain power (hereafter 
called the y "-truncated model, n signifying the highest power of Y after truncation). At 

the outset we state that the reduced model has bounded solution over the entire instability 
domain. 

In the case of the numerical solutions obtained from Eq. (|20|) , we have found that while 
truncating the ^-expansion in Eq. (23) at Y 3 is sufficient to give a bounded solution for 
the range of supercritical bifurcation (0.4713 < a < 0.512), it fails beyond this. On the 
other hand, truncating the V-expansion at Y 5 does extend the domain of bounded solution 
(0.4538 < a < 0.512) which matches with the domain where the quintic TDGL equation 
works. 

Figure 5 shows the plots of the solutions to the secular equations (Eqs. B5 and B6), the 
reduced model (Eqs. 17-18) and the F 3 -truncated model (Eq. ^(J) for a = 0.5 and e = 10~ 4 . 
All of them are numerically almost indistinguishable. We have verified that the solutions 
obtained by these three methods show that the amplitude of the limit cycle scales with ea. 

Figure 6 shows the solutions obtained by the three ways for a = 0.468 and e = 10~ 4 
where, this time, the truncation is at Y 5 . The agreement between the three solutions is 
reasonably good. Part of the discrepancy between the secular motion and the two numerical 
solutions can be attributed to the range of values of the various terms (fc's, ts, m's and n's 
of Appendix C) contributing to rj and v. (They range from 10~ 16 to 10 9 ). We have also 
verified that the amplitude of the limit cycle has a finite jump for this case. 

For a < 0.4538, where the quintic TDGL equation is inadequate, y 5 -truncated model 
is also inadequate. However, y 7 -truncated model has a bounded solution for the region 
0.4385 < a < 0.5120. In particular, Fig. 7 shows the plots of the solutions obtained by the 
y 7 -truncated model and that by the reduced model for a = 0.44 and e = 10~ 4 . This implies 
that seventh order TDGL equation need to be considered for this range of a. We have verified 
that Y 9 term extends the domain of bounded solutions upto a = 0.4246. From this trend, 
it appears that very high powers in Y need to be retained to cover the entire domain upto 
a = 0.333. 

The region a > 0.5120 is even more interesting. While the quintic TDGL equation is 
inadequate to give a bounded solution, we found that retaining even upto the ninth power 
in Y does not give a bounded solution for the entire range of values 0.512 < a < 0.707 (even 
for those close to a = 0.512). This suggests that this region is one where the nonlinearity is 
strong. On the other hand, as pointed out earlier the reduced model has limit cycle solution 
in this range as well. In particular, Fig. 8 shows such a plot for a = 0.63 and e = 10~ 4 . 

Thus, the results of the TDGL equation, the y n -truncated model and the reduced model 
are consistent with each other. Further, the numerical solutions via the power series in Y 
throws light on the nature of nonlinearities governing the solution in different regions of 
instability. 

5 Summary and Discussion 

We have carried out the reductive perturbative approach to the problem of steps on creep 
curve and shown that the dynamics of the system is described by a TDGL equation for 
the amplitude of the slow modes (complex order parameter) \l/ in the neighbourhood of 
the bifurcation point. Since the above derivation is valid only in the neighbourhood of the 
critical value, the expression for \& is valid only for small e. This has been exemplified by the 
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reasonable agreement between the secular motion and the numerical solution of the reduced 
model within the domain of the validity of the quintic TDGL equation. 

We have shown that both subcritical as well as supercritical bifurcations are seen in the 
instability range of the parameter a. While the subcritical bifurcation ( abrupt or 'first order' 
transition) seen over most of the values of a is consistent with our earlier calculation [Anan- 
thakrishna & Sahoo, 1981b; Valsakumar & Ananthakrishna, 1983], we see a supercritical 
bifurcation ( continuous or 'second order' transition) over a narrow mid-range of a. The 
existence of the supercritical bifurcation was missed by us earlier. The existence of both 
supercritical and subcritical bifurcations is consistent with experiments. In addition, the de- 
pendence of the amplitude and the periodicity of the strain jumps on stress and temperature 
are consistent with experiments. From this point of view, even though the derived TDGL 
equation is valid in a limited domain, the present calculation allows a direct mapping of 
these quantities. This calculation also helps us to demonstrate the complicated dependence 
of the slow modes on the original modes. This will serve as a warning to those using hand 
waiving arguments for declaring certain modes as fast modes and others as slow modes in 
modelling of such of these problems. 

Looked at from the point of view of the properties of the model, there are some interesting 
features. For the region outside the validity of the quintic TDGL equation, on the basis of a 
careful numerical analysis, we found that high order nonlinearities govern the nature of the 
subcritical bifurcation. It must be emphasized that this feature is unusual and this insight 
would not have been possible but for the numerical solutions obtained by keeping successive 
higher powers in Y, which in itself is the basis for the derivation of the TDGL equation. 

This unusual feature of the model interpreted in the language of phase transitions give 
better insight and could also find applicability. For conventional models, the 'free energy' 
is described by an expansion in power series of the order parameter. Even in the case of 
dynamic transition such as the present one, it is possible to associate a 'free energy' like 
function such that ^ = — where is the 'free energy' like function. While the 

'free energy' for 'second order' phase transitions is described by retaining up to fourth power 
in the order parameter, up to sixth power is conventionally required for 'first order' phase 
transitions (with the appropriate signs for the coefficients in the expansion). In the present 
case, however, we need to go to as high as twelfth power (or more) of the order parameter 
to cover the interval 0.3333 < a < 0.4246 (or equivalently twelfth power or higher in Y). 
In this case, including successively higher powers in Y increases the domain of description. 
On the other hand, in the interval of 0.512 < a < 0.707, even retaining upto tenth power 
in Y (which is the highest power we have checked) does not work. This must be contrasted 
with the existence of bounded solution for the reduced model over the entire interval. This 
feature is rather unusual and, to the best of our knowledge, we are not aware of any other 
model which exhibits this property. 
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Appendix A 

In this appendix we give the expressions for the coefficients appearing in the TDGL equation 
up to the quintic term, i.e. the expressions for fi, r\ and v in the equation J^\I/ = e/i^ + 
Tyj^l 2 ^ + z/|\I>| 4 \I>. They are given by 



and 
where 

and 



r] = 2P 2 (k 1 + k 2 + k* 2 + fc*) + 3P 3 , 
v = 2P 2 Q 2 + 3P 3 Q 3 + 4P 4 Q 4 + IOP5, 



2w 



Q 2 = (h + fc^)4 + (k 2 + fc^ + (A;* + k 3 ){e 2 + t A ) + m 5 + m 6 + m* + m* 7 , 

Q 3 = 2(h + k* 3 )(kl + k 2 + k\ + fc 3 ) + (k 2 + k 2 ) 2 + £ 2 + £ 3 + 2t 3 + t 4 , 
Q 4 = kl + k 3 + 3{h + k 2 + k* 2 + k* 3 ). 
The expressions for £ n 's and Cn's used in Eq.(A5) are given by 



6 = 



box 



a 3 
b X 



a 
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2 _ A 2 



2a 2 — a(a — b x) , 
- v« - box) (a A - 2« 2 ), 



16a 7 



6 = 



a 2 - A 2 



16a 9 



and 



c 2 = 



a(a 2 - 5A 2 ) + 8Aa 2 
aA(3a 2 - 7A 2 ) - 2a 2 (a 2 - A 2 ) 
c 



46 a 3 



(a - A ), 



A 

C3 = 2 C2, 



a^ 



C5 



. a 2 -5A 2 
C4 = 4a^ 2 ' 
A(3a 2 - 7A 2 ) 



4a 6 



C2, 



where 



A = a - 2b oX - 

Expressions for /c's, ts, m's and n's are given in Appendix C. 



(Al) 
(A2) 

(A3) 
(A4) 

(A5) 

(A6) 
(A7) 

(A8) 

(A9) 
(A10) 
(All) 

(A12) 
(A13) 
(A14) 
(A15) 

(A16) 



16 



Appendix B 

In this appendix the equations describing the secular motion are derived for the case when 
the order parameter, \I>, is found from the quintic TDGL equation (Eq. 44). 
From Eq. (39), 

Y(t) = ab (^e tuJT + c.c), (Bl) 

and 

Z(t) = (-c + iu)^e tuJT + c.c, (B2) 

where the complex amplitude is \& = e^fa + efa + ••• + €^fa with fa, fa given in terms of 
functions of fa and phase factors (see Appendix C). On the other hand, to O(e^), 

V = e^fa, (B3) 

which, upon inverting, gives fa in terms of 

fa = e-3|*|e rflT . (B4) 

Substituting the fa, V's-expressions from Appendix C in Eqs. (Bl) and (B2) and, then, 
using Eq. (B4) leads to obtaining the secular equations given by: 

Y = 2ab P + , (B5) 

and 

Z = -2c P+ -2uP-. (B6) 

where 

P± = mF±(n c T) + K±\y\ 2 F ± (2n c T + 6±J + \y\ 2 G±(k 2 ) 

+eLi |*| F ± (SI c t- 9 h ) + L± \^\ 3 F ± (3fi c r + 0±J + L 3 |*| 3 F ± (SI c t - 6 h ) 
+eM± |*| 2 F± (2fi c r + 0±J + e |*| 2 G ± (m 2 ) + M± |^| 4 F ± (4fi c r + ^ 4g ) 
+M± |*| 4 F ± (2fi c r + + |^| 4 G ± (m 6 ) + e 2 JV X |*| F± (Q c r - 9 ni ) 

+eN& |*| 3 F ± (3Q c t + 6f l2A ) +eN 3 \^ F ± (Q c r - n3 ) + N± F ± (5Q c r + 9± g ) 
+N± I vl/ 1 5 F ± (3fi c r + 0± 68 ) + N± |*| 5 F± (n c r - 0± 7 ) , (B7) 

with 



F+(0) = 


cos(4>), 


(B8) 


F_(0) = 


sin((f>), 


(B9) 


G+(u) = 


Real(u), 


(BIO) 


G_(«) = 


Im(u), 


(Bll) 




uj + Q, 


(B12) 


v± = 


\Vi±Vj\, 


(B13) 


Vi = 


\ v *l 


(B14) 


9$. = 


/ Real(vi ± 


(B15) 



where V stands for either K, L, M or N, while v stands for either k, £, m, or n correspond- 
ingly. Expressions for K, L, M and for N can be obtained from fc's, €s, m's and n's which 
are given in Appendix C. 
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Appendix C 

Expressions for ip 2 , tys, V^j ^5 (with a factor e lulT ) are 



2iu!T 



h i>le 2wr + k 2 \^\ 2 + k 3 ^f e~ 2iuJT , (CI) 
ti ft e" lwr + £ 2 ^ 3 e 3 ^ r + | 2 ft e~ iu}T + h i/>? e' 3iuJT , (C2) 
mi e 2i " r + m 2 |^ 2 | + m 3 itf e~ 2lWT + m 4 Vi e 4tuJT + m 5 |^| 2 V? e 2 ^ r 
+m 6 |^| 4 + m 7 iV'ilVf e" 2 ^ + m 8 e" 4 ^, (C3) 
m V? e~^ T + n 2 V? e 3lwr + n 3 Vi e~^ r + n 4 e~ 3iWr + n 5 e 



',5 „5jwr 



+n 6 |^i|V? e 3l " r + n 7 |^i| 4 ^i e" icJr + n 8 |ViNf e" 3l " r + n 9 ^ e~ 5luJT (CA) 
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7. Pi, 

6L0 
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4a; 
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n 2 = -—{2P 2 [f 1 (k 1 + k$)+m 1 + m* 3 ] + 3P 3 ei-2 f d 2 -ii%}, 



n 3 = — {2P 2 [£ 1 (k 2 



fc*) + £*(fc* + h) +1711+1712+1712+ m 3 ] 

+3P 3 (£* + 2£ 1 ) -rf l x - 2// £ 3 } , 



n 4 



— {2P 2 + fc 3 ) + m\ + m 3 \ + 3P 3 h - tf£% + ~ , 

4cj 



n 5 



= {2P 2 [m 4 + ml + (h + fc*) (4 + 3D] + 3P 3 [(*i + fc 3 ) 2 + £ 2 + i 

+4P 4 (A; 1 + fc*) + P 5 } , 

n 6 = -^{2P 2 [m 4 + m 5 + m^ + ml + £l(k 1 + k*) + (e 2 +t 4 )(k2 + k*)} 
2ui 
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i\ + 2{l 2 + t A ) + (h + klY + 2(h + fc*)(fc 2 + k* 2 



+4P 4 [k 2 + k* + 3(h + fc*)] + 5P 5 -3^ 2 } , 
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m 6 + m 7 + {t\ + h)(h + fc*) + £ 3 (fc 2 + fc*) + £*(fc* + fc 3 )] 



(fc 2 + k* 2 Y + 2(fc* + fc 3 )(fci + k 2 + k* 2 + kl) +t 2 + t 3 + 2t 3 + £ 4 



+4P 4 [fc! + fc* + 3(fc* + fc* + fc 2 + k 3 )\ + 10P 5 -(rj + 2r}*)£ 3 } , 
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+3P 3 [(fc* + fc 3 ) 2 + 2(fc 2 + fc*)(fc* + k 3 ) + 2[l\ + £4) + h 
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6u 



n 9 = —{2P 2 [ml + m 8 + (t2 + £4)(K + k 3 }+3P 3 [(k* + k 3 ) 2 + e 4 + £ 4 
+4P 4 (fc* + fc 3 ) + P 5 } . 

Note that expressions for P 2 ,...,P 5 are given in Eq. (A5) in Appendix A. 
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Figure Captions 



Fig. 1 Plot of — versus a with b = 10 4 . The hatched portion shows the instability region 

Cq 

which extends between a — | and a = ^ approximately. 

Fig. 2 Plots of 771 and v\ versus a (b = 10~ 4 ). Note that rji is negative (but finite) in 
the interval 0.4713 < a < 0.5210. Similarly, v\ is negative (but finite) in the interval 
0.4538 < a < 0.5097. 

Fig. 3 Plots of and P versus a in the region of supercritical bifurcation. 

Fig. 4 Plots of 1^1 2 and P versus a in a subregion of subcritical bifurcation, i. e. the region 
where quintic TDGL equation holds. 

Fig. 5 Plots of the (numerical) limit cycle solutions obtained from the secular equations 
(dotted line), the reduced model (dashed line) and the F 3 -truncated model (crossed 
line) for a = 0.5 and e = 10~ 4 . 

Fig. 6 Plots of the (numerical) limit cycle solutions obtained from the secular equations 
(dotted line), the reduced model (dashed line) and the F 5 -truncated model (crossed 
line) for a = 0.468 and e = 10" 4 . 

Fig. 7 Plots of the (numerical) limit cycle solutions obtained from the reduced model 
(dashed line) and the F 7 -truncated model (crossed line) for a = 0.44 and e = 10~ 4 . 

Fig. 8 Plot of the (numerical) limit cycle solution obtained from the reduced model for 
a = 0.63 and e = 10" 4 . 
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